↜ Back to index Introduction to Numerical Analysis 1

Part b—Lecture 4

We discuss another important partial differential equation, the wave equation.

The wave equation

This time we consider a string stretched along the x-axis and attached at points 0 and 1. At time t = 0 and point x ∈ [0,1], the string position above the x-axis is u_0(x) and the string is moving in the vertical direction with velocity v_0(x). We are interested in the position of the string at points x \in (0,1) at a later time t > 0, called u(x,t). The velocity of the string in the vertical direction is then given by the time derivative u_t(x, t).

Assuming that the string is under a constant tension force, Newton’s law implies that u is a solution of the wave equation:

u_{tt}(x, t) = u_{xx}(x, t), \qquad 0 < x < 1, t > 0,\\

This time we have u_{tt}, the second partial derivative of u with respect to t. For u to be a solution of the equation, it needs to be twice differentiable both with respect to x and t and the equality must hold at all x \in (0, 1) and t > 0.

Example 1. Let \omega \in {\mathbb{R}} be fixed. We check that the function u(x,t) = \cos(\omega t) \sin (\omega x) is a solution of the wave equation.

We first compute the partial derivatives using the chain rule: \begin{aligned} u_t(x,t) &= -\omega \sin(\omega t) \sin (\omega x),\\ u_{tt}(x,t) &= -\omega^2 \cos(\omega t) \sin (\omega x),\\ u_x(x,t) &= \omega \cos(\omega t) \cos (\omega x),\\ u_{xx}(x,t) &= -\omega^2 \cos(\omega t) \sin (\omega x). \end{aligned} We see that u_{tt} = u_{xx} everywhere and so this function is a solution of the wave equation.

Unlike the heat equation, the solution does not become small as time increases. Instead, the solution oscillates with angular frequency \omega.

Example 2. The wave equation is also a linear partial differential equation:

If u and v are two solutions, then u + v is also a solution.
If u is a solution, then C u is a solution for any constant C \in {\mathbb{R}}.

Exercise 1. Here are a few examples of solutions of the wave equations. Check that they satisfy the wave equation for all x and t.

u = C for any constant C.
u = C_1x + C_2t for any constants C_1, C_2.
u = \cos(\omega t) \cos (\omega x) for any \omega \in {\mathbb{R}}.
u = \sin(\omega t) \cos (\omega x) for any \omega \in {\mathbb{R}}.

As in the case of the heat equation, we need to prescribe extra conditions to have a unique solutions.

We need to prescribe:

initial condition: u(x, 0) = u_0(x), the position of the string at the initial time t = 0, and u_t(x, 0) = v_0(x), the velocity in the vertical direction at t = 0. Here u_0 and v_0 are given functions.
boundary condition: u(0, t) = a, u(1, t) = b, where a, b are given constants, prescribe how far above the boundary points is the string fixed. This is again called the Dirichlet boundary condition.

The initial value problem for the wave equation with Dirichlet boundary condition:

\left\{ \begin{aligned} u_{tt}(x, t) &= u_{xx}(x, t), && 0 < x < 1, t > 0,\\ u(x, 0) &= u_0(x), && 0 < x < 1, \\ u_t(x, 0) &= v_0(x), && 0 < x < 1, \\ u(0, t) &= a, && 0 < t, \\ u(1, t) &= b, && 0 < t. \end{aligned} \right.

Numerical method

The setup is the same as with the heat equation. We subdivide the domain (0, 1) into M \geq 1 subintervals (x_{k-1}, x_k) with x_k = k h, h = 1/M. We also choose a time step \tau and define the times t_n = n \tau. We look for u_{n, k} that approximates u(x_k, t_n), k = 0, \ldots, M, n \geq 0.

Since both sides of the wave equation have the second derivative, we use the same idea of the central difference and get a scheme \frac{u_{n+1, k} - 2u_{n, k} + u_{n-1, k}}{\tau^2} = \frac{u_{n, k - 1} - 2 u_{n, k} + u_{n, k + 1}}{h^2}, \qquad n = 0, 1, \ldots, \quad k = 1, \ldots, M-1.

We can express u_{n+1, k} in terms of u_{n, k - 1}, u_{n, k} and u_{n, k + 1} and u_{n-1, k}. Therefore we can advance the solution from timestep n to timestep n+1 if we also know the solution at timestep n-1.

The values u_{0, k} for k = 1, \ldots, M-1 are given by the initial data u_0(x_k), and the values u_{n, 0}, u_{n, M} for n =0, 1, \ldots are given by the boundary data a and b.

Since we need to know the solution at 2 previous time steps, we have to take extra care to find the values at step n = 1 since we do not know the values at step n = -1. Here we use a slightly different finite difference, using the initial velocity v_0:

u_{1, k} = u_0(x_k) + \tau v_0(x_k) + \frac{\tau^2}{2h^2} \big(u_0(x_{k-1}) - 2 u_0(x_k) + u_0(x_{k+1})\big), \qquad k = 1, \ldots, M - 1.

This is motivated by the fact that the solution of y''(t) = c with constant c is y(\tau) = y(0) + \tau y'(0) + \frac{\tau^2}2 c.

Algorithm.

Set the initial values at n = 0: u_{0, k} = u_0(x_k), \qquad k = 0, \ldots, M.
Compute the values at n = 1: \begin{aligned} u_{1, 0} &= a,\\ u_{1, k} &= u_{0,k} + \tau v_0(x_k) + \frac{\tau^2}{2h^2} \big(u_{0,k-1} - 2 u_{0,k} + u_{0,k+1}\big),\\ &\qquad\qquad\qquad k = 1, \dots, M-1\\ u_{1, M} &= b. \end{aligned}
For every n = 1, 2, \dots, find the “new” values: \begin{aligned} u_{n+1, 0} &= a,\\ u_{n+1, k} &= 2u_{n, k} - u_{n-1, k} + \frac{\tau^2}{h^2} (u_{n, k - 1} - 2 u_{n, k} + u_{n, k + 1}), && k = 1, \ldots, M-1,\\ u_{n+1, M} &= b. \end{aligned}